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I. INTRODUCTION 

Nucleon electromagnetic form factors are fundamental quantities and reveal important information on the spatial 
distribution of charge and magnetization within a nucleon [IHS] . Understanding the nucleon electromagnetic structure 
in terms of the underlying quark and gluon degrees of freedom of quantum chromodynamics is a challenging task 
which has attracted the attention of both theory and experiment for many years. 

For a long time, the overall trend of the experimental results for small and moderate values of the momentum 
transfer Q 2 = —q 2 could be described reasonably well by phenomcnological (dipole) fits 



G p E (Q 2 ) 



GIAQ 2 ) giaq 2 ) 



~ (l + Q 2 /m 2 D )- 2 , 
G n E (Q 2 )~0, (1) 
with mo 0.84 GeV and the magnetic moments 

~ 2.79 , fi n ~ -1.91 , (2) 

in units of nuclear magnetons. More recently, the improved accuracy of the experimental data allows us to see clear 
deviations from this dipole behavior in the region of low and intermediate Q 2 . This has led to a significant amount 
of theoretical work aimed at describing these form factors, such as dispersion theory analysis [ME], vector meson 
exchange/dominance [9], [10] and polynomial fits based on the Kelly parametrization [11] . We will apply a simplified 
variant of the latter to lattice results in this paper. Interest in these form factors has been revived over the last 10 
years by experiments at Jefferson Laboratory which found an unexpected dependence of the nucleon's electric and 
magnetic form factors on the momentum transferred to the target nucleon [12H14] , More recently, a measurement 
of the Lamb shift in muonic hydrogen |15j has produced a result for the electric radius of the proton that is several 
sigma below the PDG (CODATA) value [THIE]. At the same time, a new high- precision determination of the proton 
form factors from ep-scattering experiments at MAMI has been reported [IB], which confirms the traditional results 
for the electric (and magnetic) mean square radius. 

From a lattice perspective, it is common to evaluate fundamental observables like the charge radii and anomalous 
magnetic moments in order to make a comparison with experimental and theoretical results. A feature of any 
lattice simulation is that the quark mass is an input parameter, hence it is possible to map out the form factors 
not only as a function of Q 2 , but also m\. Baryon charge radii and magnetic moments are of particular interest 
in this case as predictions from chiral perturbation theory (ChPT) indicate that these quantities should provide an 
excellent opportunity to observe the chiral non-analytic behavior of QCD [TM2"1] . An additional advantage of a lattice 
simulation of nucleon electromagnetic form factors is that since they are performed at the quark level, it is possible to 
determine the individual up and down quark contributions, providing valuable insights into the distribution of charge 
and magnetization within a nucleon. 

These issues are now beginning to be addressed in modern lattice simulations [221 [2"51 - f2"9"] (see also [3D] for a review) . 
A common feature of present lattice simulations with unphysical quark masses is that they tend to underestimate 
the experimental and phenomenological results for the radii and magnetic moments of the nucleon. As mentioned 
above, predictions from ChPT indicate that these observables should exhibit a dramatic non-analytic dependence on 
the quark mass close to the chiral limit, however such features have yet to be seen clearly in a lattice simulation with 
dynamical quarks. 

In this paper, we will confront these issues through simulations with pion masses as low as ~ 180 MeV. In our 
analysis, we place a strong emphasis on addressing the systematic errors present in a lattice simulation, such as finite 
volume and lattice spacing effects. We also consider the effects of finite momentum resolution in lattice determinations 
of form factors through the use of several parametrizations of the momentum dependence. 

II. LATTICE SETUP AND METHODS 

Below we briefly describe our lattice setup and the methods that we have used to compute the nucleon form factors. 

A. Simulation parameters 



We perform our simulations with two flavors of non-perturbatively 0(et)-improved Wilson (Clover) fermions and 
Wilson glue. Using these actions, we have generated gauge field configurations with the parameters given in Table [TJ 
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p 


# 




N 3 


x T 


m„ [GeV] 


a [fm] 


L [fm] 


ATtraj 


5.20 


1 


0.13420 


16 3 


x 32 


1.40 


0.083 


1.3 


0(5000) 


5.20 


2 


0.13500 


16 3 


x 32 


0.99 




1.3 


0(8000) 


5.20 


3 


0.13550 


16 3 


x 32 


0.69 




1.3 


0(8000) 


5.25 


4 


0.13460 


16 3 


x 32 


1.29 


0.076 


1.2 


0(6000) 


5.25 


5 


0.13520 


16 3 


x 32 


1.00 




1.2 


0(8000) 


5.25 


6 


0.13575 


24 3 


x 48 


0.67 




1.8 


0(6000) 


5.25 


7 


0.13600 


24 3 


x 48 


0.48 




1.8 


0(5000) 


5.29 


8 


0.13400 


16 3 


x 32 


1.59 


0.072 


1.1 


0(4000) 


5.29 


9 


0.13500 


16 3 


x 32 


1.16 




1.1 


0(5500) 


5.29 


10 


0.13550 


12 3 


x 32 


0.99 




0.9 


0(4500) 


5.29 


11 


0.13550 


16 3 


x 32 


0.92 




1.1 


0(5000) 


5.29 


12 


0.13550 


24 3 


x 48 


0.90 




1.7 


0(2000) 


5.29 


13 


0.13590 


12 3 


x 32 


0.93 




0.9 


0(5500) 


5.29 


14 


0.13590 


16 3 


x 32 


0.69 




1.1 


0(7000) 


5.29 


15 
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24 3 
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0.66 
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0(6000) 


5.29 


15 
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24 3 
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0.43 
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5.29 


17 


0.13632 


24 3 


x 48 


0.31 




1.7 


0(7000) 


5.29 


18 


0.13632 


32 3 


x 64 


0.30 




2.3 


0(2700) 


5.29 


19 


0.13632 


40 3 


x 64 


0.29 




2.9 


0(2000) 


5.29 


20 


0.13640 


40 3 


x 64 


0.18 




2.9 


0(1000) 


5.40 


21 


0.13500 


24 3 


x 48 


1.32 


0.060 


1.4 


0(3500) 


5.40 


22 


0.13560 


24 3 


x 48 


1.02 




1.4 


0(3500) 


5.40 


23 


0.13610 


24 3 


x 48 


0.72 




1.4 


0(4000) 


5.40 


24 


0.13625 


24 3 


x 48 


0.62 




1.4 


0(6000) 


5.40 


25 


0.13640 


24 3 


x 48 


0.50 




1.4 


0(2500) 


5.40 


26 


0.13640 


32 3 


x 64 


0.49 




1.9 


0(2500) 


5.40 


27 


0.13660 


32 3 


x 64 


0.28 




1.9 


0(2800) 


5.40 


28 


0.13660 


48 3 


x 64 


0.26 




2.9 


0(2200) 



TABLE I. Overview of our simulation parameters where we have used the Sommer parameter with ro = 0.5 fm to set the 
physical scale. 



where we have used the Sommer parameter with ro = 0.5 fm to set the physical scale |31j . Summarizing these 
parameters, we see that our four values of /3 = 5.20, 5.25,5.29,5.40, correspond to lattice spacings in the range 
0.06 < a < 0.1 fm, allowing for the approach to the continuum limit to be assessed, while a range of lattice volumes 
(0.9 < L < 3.0 fm) enable us to study finite size effects in our simulations. Finally, our pion masses now reach well into 
the chiral regime, down to m, ~ 180 MeV, allowing us to investigate the applicability of different ChPT approaches 
around and above the physical pion mass, and to search for chiral non-analytic behavior in our results. When 
computing correlation functions on these configurations, we generally over-sample using up to 4 different locations of 
the fermion source on a single configuration. We then use binning to obtain an effective distance of 20 trajectories. 
We find that beyond this, the size of the bins has little effect on the error, which indicates residual auto-correlations 
are small. 



B. Extraction of form factors 



On the lattice, we determine the form factors F\(Q 2 ) and ^(Q 2 ) by calculating the following matrix element of 
the electromagnetic current 



(p>, s'|j*(0)b, a) = U{ P >, s') 



^ Fl ( q 2 ) + ta^^F 2 ( q 2 ) 



U(p, s) , 



(3) 



where U(p, s) is a Dirac spinor with momentum p and spin polarization s, q = p' — p is the momentum transfer, 

is the nucleon mass and j M is the electromagnetic current. The Dirac, F±, and Pauli, i*2, form factors of the proton 
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are obtained by using 



2 1 _ 

$ = 3«7m m - 3^7^ ( 4 ) 



between proton states. The isovector form factors are also obtained from proton states, but with the current 

Fn = «7/*« - <*7/*d • ( 5 ) 

Similarly, we used the isoscalar current uj^u + dj^d for the computation of isoscalar form factors. 

In electron scattering, it is common to rewrite the form factors F% and F 2 in terms of the electric and magnetic 
Sachs form factors, 

G E (Q 2 )=F 1 (Q 2 )- J -^F 2 (Q 2 ), 

G M (Q 2 )^F 1 (Q 2 ) + F 2 (Q 2 ), (6) 

as then the (unpolarized) cross section becomes a linear combination of squares of the form factors. 

For, e.g., the proton F^\o) = G^(0) = 1 gives the electric charge, while G$(0) = m (p) = 1 + « (p) gives the 
magnetic moment, where F^ (0) = kW is the anomalous magnetic moment. For a classical point particle, both form 
factors are independent of Q 2 , so deviations from this behavior tell us something about the extended nature of the 
nucleon. 

In our lattice study, we use the standard proton interpolating field for a proton with momentum p 

B a (t,p)= J2 e-'*%k<(xH( x X C ri(hA( x )> 

X,X4=t 

B a (t,p) = j p -%kfy(x)(CK)^(x)u k a (x), (7) 



X,X4— t 



where C is the charge conjugation matrix, i , j ,k are color indices and a , /3 , 7 are Dirac indices. 

In order to improve the overlap of these interpolating fields with the ground state proton, we employ two improve- 
ments: Jacobi smearing and non-relativistic projection. The latter of these has the additional advantage that we only 
need to perform 2x3 inversions rather than the usual 4x3, since we only consider the first two Dirac components. 

The matrix elements in Eq. ^ are obtained from ratios of three-point to two-point functions, 



R(t,T;p',p;0) = 



C% t (t,T;p',p,Q) 
C2 P t(t,p') 



C 2p t{T,p')C2 p t{t,p')C 2p t{t -T,p) 
C 2 pt(T,p)C 2 pt(t,p)C 2pt (t - T,p') 



(8) 



for large time separations, « t « f < ^Lt, where Lt is the temporal extent of our lattice. The nucleon two- and 
three-point functions are given, respectively, by 

C 2p t(T,p) = Tr [T unpol (B(T,p)B(0,p))] , 
Ci pt (t,r;p\p,0) = Tr [T(B(t,p')O(q,T)B(0,p))] . (9) 

Here t and r are the Euclidean times of the nucleon sink and operator insertion, respectively, p' (p) is the nucleon 
momentum at the sink (source), and O is the local vector current 

0„% t) = £ e^ s q(x, r)^q(x, r) . (10) 

X 

The trace in Eq. ^ is over spinor indices and the T matrix in the three-point function is one of 

r unpo i = ^(1 + 74) , (11) 

Ti = ^(l + 74)i757i , (12) 
F 2 = -(1 + 74)17572 ■ (13) 
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We simulate with three different sink momenta p ' 

= (0, 0, 0), (1, 0, 0), (0, 1, 0) . (14) 

Finally, we use 17 different momentum transfers q — p' — p. Equations with identical values of virtual momentum 
transfer q 2 are combined to return the optimal statistics available at each working point. This procedure is outlined 
in more detail in Ref. [25] . 

Note that quark line disconnected contributions to the three-point function in Eq. ([9]), which are relevant for the 
flavor singlet observables but cancel out in the isovector case, have not been included in our study. 

III. NUMERICAL RESULTS 

In the following, we present and discuss in some detail our numerical results for the Dirac and Pauli nucleon form 
factors. 

A. Q 2 -dependence of Fi and F2 

In Figs. [lJtoH we provide an overview of our results for the Dirac and Pauli form factor^] in the isovector (u — d) 
and isosinglet (u + d) channels, for three different ranges of m n , including all ensembles specified in Table [i] Lighter 
pion masses are indicated by darker colored points. 

For comparison, we also show in each case the parametrization of the Q 2 "dependence of the experimental data 
obtained by Alberico et al. |32] as gray error bands. This parametrization has been originally performed for Sachs 
electric and magnetic form factors for the proton and the neutron, G p ^ n M , using Kelly's parametrization ansatz |llj 
for G P E M and G^ /7 and a Galster parametrization for G^, with in total 14 parameters. Since the parameters are 
strongly correlated, we have employed the full error correlation matrices provided at |33j for the error propagation. 
For our purposes, we consider the resulting parametrization of the Q 2 -dependence of F^ d as a reasonably faithful 
representation of the experimental data, at least for not too large Q 2 , and will use it as such throughout this work. 
We note that since the quality and availability of experimental results at larger values of the momentum transfer 
above ~ 1 GeV 2 decreases, in particular for G^, the shown error bands might significantly underestimate the actual 
uncertainties in this region. If not stated otherwise we will consider the proton form factors in the following and omit 
the superscript p for notational simplicity. 

In the cases of F" in Figs. |l|and || the normalization at Q 2 = is fixed, and we clearly see that the slope of the 
lattice data is significantly smaller than that of the parametrization. It is interesting to observe, however, that the 
data points systematically move towards the physical result as the pion mass decreases. Concerning F^~ d in Fig. [2] 
it seems at first sight that the lattice data points, which show little dependence on m ff , are in rough agreement with 
experiment over a wide range of Q 2 . This can be quite misleading, as we will see in more detail in the following 
sections: Not only is the slope of the lattice data points too small, but the lattice results for F2~ d (Q 2 = 0) = n u ~d 
(obtained from extrapolations in Q 2 ) are also significantly below the experimental value. In combination, one naturally 
finds that the lattice and the experimental results do overlap in a certain range of the momentum transfer, however 
without implying a general agreement for all Q 2 . Finally, our results for in Fig. [4] turn out to be compatible 

with zero within errors for practically all accessible values of Q 2 , with the exception of a small number of data points 
at lower pion masses and low momentum transfers, in particular the lowest Q 2 ~ 0.2 GeV 2 , showing a slight trend 
towards negative values and the experimental error band. 

Before engaging in a more detailed study of the dependence of F12 on Q 2 and also on the pion mass, we briefly 
address two interesting questions that have been discussed before in the literature and that can be addressed on the 
basis of ratios of form factors. 

The first deals with potentially different Q 2 -slopes of F\ for the up and the down quarks. In coordinate/impact 
parameter space [34] [35] , such different slopes would correspond to differently shaped quark density distributions and 
therefore provide important information about the inner structure of the nucleon. Figure [3] gives an overview of our 
results for the ratio F d /F™ as a function of Q 2 , where we have included all lattice ensembles, and where darker colors 
correspond to lighter pions on the lattice. The parametrization of the experimental data is, as before, illustrated by 



1 As in our earlier study in I23| , we have normalized the results for F2 such that the anomalous magnetic moment is given in units of the 
physical nuclear magneton, e/(2m^ yB ). 
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FIG. 1. Dirac form factor Fi(Q 2 ) in the isovector channel. FIG. 2. Pauli form factor J2(Q 2 ) in the isovector channel. 

All ensembles are included, and darker colors correspond to All ensembles are included, and darker colors correspond to 

lighter pion masses. The gray shaded band represents the lighter pion masses. The gray shaded band represents the 

parametrization by Alberico et al. [32] of the experimental parametrization of Ref. [32] of the experimental data, 
data. 



the error band. While the data points at the largest pion masses show only a small dependence on Q 2 , they show a 
systematic downwards trend towards ratios Ff /F" < 0.5 for Q 2 > 0.5 GeV 2 as the pions get lighter, thereby moving 
closer to the experimental error band. In the limit Q 2 — > 0, the experimental result flattens off considerably, which 
we will discuss in more detail in section III C below on the basis of the separate mean square radii for up and down 
quarks and different parametrizations of the lattice data. 

The second question concerns the scaling of F 2 /Fi at intermediate to large <2 2 -values. Perturbative QCD suggests 
that Q 2 F 2 /Fi ~ const, as Q 2 — > oo, up to logarithmic corrections oc InQ 2 [36-38 ]• Recently, supported by new 
measurements of at JLab Hall A, it has been noted that F2/F1 approximately scales as a constant already in 
an intermediate range of Q 2 = 1.5, ... ,3.5 GeV 2 , separately for up and for down quarks [35]. We show our results 
together with the parametrization of [32] for these ratios in Fig s. [6| and [7] Clearly, the uncertainties and the scatter 
of the lattice data, in particular for the up quark case in Fig. |6J make it difficult to draw any strong conclusions. 
Nevertheless, the lattice data points in the different m^-ranges are overall compatible with a flat Q 2 -dependence above 
~ 0.5 GeV 2 , which is most clearly seen for the down quarks in Fig. [7] at the largest pion masses. While the latter 
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FIG. 3. Dirac form factor Fi(Q 2 ) in the isosingiet (u + d) 
channel. All ensembles are included, and darker colors cor- 
respond to lighter pion masses. The gray shaded band rep- 
resents the parametrization by Alberico et al. [32] of the 
experimental data. 
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FIG. 4. Pauli form factor F 2 (Q 2 ) in the isosingiet (u + d) 
channel. All ensembles are included, and darker colors 
correspond to lighter pion masses. The gray shaded band 
represents the parametrization of Ref. [32] of the experimen- 
tal data. 



are about a factor of two below the experimental band, we see a clear upwards trend as lower pion masses are being 
approached. These trends are less clear in the case of the up quarks in Fig. [6j for which the data points are, however, 
generally closer to experiment. 

In any case, more quantitative conclusions with respect to these interesting questions will have to be based on 
precise lattice data at low pion masses that extends up to and beyond squared momentum transfers of Q 2 ~ 2 GeV 2 . 

To conclude this section, we note that a study of potential systematic uncertainties, as well as the pion mass 
dependence of the form factors at fixed Q 2 , is given further below in section jVj In short, we do not see any significant, 
systematic effects due to contributions from excited states (section VA|, or the finite lattice spacing (section [VB ), at 
least within statistical uncertainties. Although the lattice data points show an approximately linear dependence on 
or m 2 at fixed Q 2 , we find that simple linear extrapolations to the physical point would not lead to an agreement 
with the results from experiment and phenomenology (section VB). We therefore conclude that a non-trivial pion 
mass dependence has to set in between the lowest accessible lattice pion masses of ~ 180, . . . , 260 MeV and mP hys . 
This will also be studied in greater detail on the basis of the mean square radii and anomalous magnetic moments in 




0.1 - 
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G 2 [GeV 2 ] 

FIG. 5. The ratio Ff / F™ of down to up quark contributions to the Dirac form factor. All ensembles are included. The darker 
colors correspond to lighter pion masses. The gray shaded band represents the parametrization of Ref. [35] of the experimental 
data. 



section IIVI 



B. Parametrizations of the Q 2 -dependence 

We now turn to analytical parametrizations of the Q 2 -dependence. These will not only allow us to interpolate 
between the discrete values of Q 2 , but in particular to extrapolate our results for F2 to the forward limit in order to 
extract the anomalous magnetic moment. Furthermore, well-chosen parametrizations are important to obtain more 
realistic estimates for the mean square radii from the slopes of the form factors at Q 2 = 0, 



2 6 dF t (Q 2 



Fi(0) dQ 2 



(15) 

5 2 =0 



However, we note that the parametrizations unavoidably introduce some model dependence into the analysis. 

In the following, we compare different ansatze for the Q 2 -dependence of the form factors. A common ansatz for the 
Dirac form factor is a dipole, 

F '& ) = v/wUr (16) 

with, e.g., F±' u ~ (0) = 1, where the dipole mass mu is a free fit parameter. The corresponding mean square radius 
is then given by the squared inverse dipole mass, (r 2 )i = 12/m 2 D . A more flexible parametrization is obtained with a 
more general polynomial in the denominator, 

F 1 (Q 2 )= ^ F 'i°l _ 4 , (17) 

1 + ci 2 q z + C14Q 4 

with C12 and C14 as free fit parameters. Here the mean square radius is obtained from (r 2 ) x = 6c 12 . The latter form 
was already employed in Ref. [40] . and it also allows for a matching to a simple vector meson exchange ansatz, as 
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FIG. 6. The form factor ratio F2/(k,F\) for up quarks. All en- 
sembles are included. The darker colors correspond to lighter 
pion masses. The lattice data points have been obtained us- 
ing the experimental values for k„ in the ratio. The gray 
shaded band represents the parametrization of Ref. [32] of 
the experimental data. 
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FIG. 7. The form factor ratio F2/(kFi) for down quarks. 
All ensembles are included. The darker colors correspond 
to lighter pion masses. The lattice data points have been 
obtained using the experimental values for hid in the ratio. 
The gray shaded band represents the parametrization of 
Ref. [32] of the experimental data. 



will be discussed below in section [HI D[ Similarly, for the Pauli form factor F 2 , one could employ a simple dipole or 
tripole form 



F 2 (Q 2 ) = 



F 2 (0) 



(l + Q 2 /m 2 p )P 



(18) 



where p = 2 or 3, and F 2 (0) and the pole mass m p are the fit parameters. In this case the Pauli radius is given by 
(r 2 )2 = 6p/rrip. Alternatively, a more general polynomial in the denominator leads to a three-parameter ansatz of the 
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FIG. 8. Parametrization of the Q -dependence of the isovector Dirac form factor lattice data for selected ensembles. The 
shaded bands represent the 2-parameter fits to the lattice data points based o n Eg . (17 1. We also show the corresponding 
values of \ 2 /DOF. For comparison, the 1-parameter dipole fits based on Eq. (16 1 are indicated by the dashed lines, with 
\ 2 /DOF given in parentheses. 



form [30] 



F 2 (Q 2 ) = 



F 2 (0) 



1 + c 22 Q 2 + c 26 Q 6 



(19) 



for which the mean s qua re ra dius is (r 2 ) 2 = 6c 22 . 



denominators of Eqs. (17) and |l9| ensure that Fi(Q' 2 ) 



We note that the choices for the highest powers of Q 2 in the 



— — >~ 1/(Q 2 ) 1+1 , as expected from perturbative QCD [36] , 

Typical results for parametrizations of F^~ d {Q 2 ) for selected ensembles are displayed in Fig. [8] We observe that 
according to the \ 2 /DOF, the less restrictive polynomial ansatz seems to describe the data significantly better. Apart 
from that, the main difference between the fits based on Eqs. (16) and p"7| is the smaller slope and the broadening 
of the error bands in the case of the 2-parameter polynomial ansatz in the region of larger Q 2 where no data points 
are available. 

In the case of F^~ d , we are comparing a tripole (p — 3) ansatz, Eq. (18), with the more general polynomial 
parametrization, Eq. (19), in Fig. [9] With respect to the x 2 /DOF, the 3-parameter form doesn't have any advantage. 
However, due to the somewhat stronger broadening of the error bands at lower and larger values of Q 2 in regions where 
no data points are available, we consider the polynomial ansatz in general to be less biased and hence to provide more 
realistic uncertainties with respect to extrapolations, in particular to Q 2 — 0. A more quantitative comparison of 
different parametrizations of ^(Q 2 ) will have to be based on precise data over a broader range of Q 2 , e.g., employing 
partially twisted boundary conditions |41L I42j to access very small momentum transfers. 

In the following sections, we will argue further on the basis of a matching to a simplistic vector meson exchange 
ansatz, and the extracted Dirac and Pauli radii, that the polynomial ansatze in Eqs. (17) and (19) provide a more 
consistent description of the data. We will therefore consider them in the following as preferred compared to the 
standard dipole and tripole forms in Eqs. (16) and (fl8|). 



A collection of numerical results for the mean square radii and anomalous magnetic moments, obtained from the 
polynomial parametrizations, is provided in Appendix [XJ for all ensembles listed in Table [TJ 
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shaded bands represent the 3-parameter fits to the lattice data points based on E g. (19 1. We also show the corresponding 
values of x 2 /DOF. For comparison, the 2-parameter tripole fits based on Eq. (18 1 (p = 3) are indicated by the dashed lines, 
with x 2 /DOF given in parentheses. 



Dirac radius 



An overview of our results for the isovector Dirac radius, 
obtained from the polynomial ansatz in Eq. 



„2\u— d 



(171. 



as a function of is provided in Fig. 
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as 



Although the results from the polynomial ansatz are somewhat larger at lower values of m n relative to the standard 
dipole fits, we find that even at the lowest accessible pion masses of 200 — 300 MeV, the lattice data points are about 
50% below the phenomenological and experimental results. It is interesting to point out, however, that we observe an 
upwards trend for < 400 MeV that does not seem to follow the otherwise rather linear pion mass dependence of the 
data points. In passing, we also note that there is a significant, so far unresolved difference between the values at the 
physical point obtained from the recent muonic hydrogen measurements |15j and the PDG |16j . Incidentally, the result 
from an earlier dispersion relation analysis of experimental form factor data [H] agrees well with the recent muonic 
hydrogen study. The observation that lattice calculations at unphysically large pion masses give mean square radii 
that are significantly below experiment has been made already in a number of previous publications, e.g. [231 I25lf2"5] . 
In combination with a detailed study of potential discretization and finite volume effects, as well as contaminations 
from excited states in section [Vj we will come to the conclusion that indeed a strong pion mass dependence has to 
set in between the physical pion mass and ~ 200 MeV. This is also in agreement with general predictions from 



chiral perturbation theory, as we will discuss below in section IV 

Apart from the normalization, the most significant difference between results from the polynomial and the dipole 
fits are the relative positions of the lattice data points for small values of m v x L < 3.4, which are most likely affected 
by finite volume effects. They turn out to be residing above the data points for larger x L in the case of the dipole 
parametrization, and below for the polynomial ansatz in Fig. 10 Since one generically expects the radius of a hadron 



to decrease as the volume decreases, we find again that the polynomial fit provides a more physical parametrization 
of our data. 

Corresponding results for the isosinglet Dirac radius are shown in Fig. 
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Overall, the data points feature small 
statistical uncertainties and show only little scatter over the full range of pion masses. Below m n < 700 MeV, we 
even find a remarkable upwards tendency, although the lattice results at m w ~ 250 MeV are still ~ 25% below the 



expected range of values of (r' 



\u+d 



0.60, 



,0.62 fm at the physical point. 
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FIG. 11. Isosinglet Dir ac r adius versus rrw, as obtained from 
fits to F 1 " +d using Eq. (fir I . 
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FIG. 12. Comparison of Dirac radii for up and down quarks FIG. 13. Ratio of the Dirac radius of down to up quarks in 
in the proton as a function of the pion mass. the proton as a function of the pion mass. 



A more detailed discussion of the pion mass dependence of the isovector and isosinglet Dirac radii will be given 
below in section |IV A| Although we cannot ex clude the presence of some discretization and finite volume effects in 



2 )i~ , our corresponding analysis in section VC does not provide any indication that they are larger than the 
present statistical uncertainties. 

In conjunction with the comparison of the slopes of F\(Q 2 ) for up and for down quarks above (see the ratio Ff/Fj 1 
in Fig. [5]), we show in Fig. 12 the corresponding Dirac radii as functions of m n . They were obtained from separate 
parametrizations of F"(Q 2 jand Ff(Q 2 ) using Eq. ( fi~7| ). We find that the mean square radii of the down quarks 
are systematically larger than those of the up quarks. The corresponding ratio in Fig. [13] is rather flat over the full 
range of pion masses, with an average of (r 2 )f/(r 2 )i ~ 1.24 ± 0.13. The observed hierarchy is in agreement with the 
experimental and phenomenological results, although the latter show a much smaller difference between the up and 
the down quark radius of just w 2 — 4%. It will be interesting to study the origin of this feature in more detail in the 
future. We also note that the substantial difference between the values obtained from the analysis of Belushkin et al. 
in [5] and the form factor parametrization of Alberico et al. [32] in Fig. [5] is, as before, at least to some extent related 
to a corresponding difference in the proton charge radius, (r 2 ) p E . 



D. Matching to a vector meson exchange ansatz 



Specific contributions from vector meson exchange, the two-pion continuum, perturbative QCD etc. to nucleon form 
factors have, for example, been investigated in [5] in the framework of a dispersion relation study of experimental 
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FIG. 14. Ratio of the lowest pole mass extracted from C12 FIG. 15. The upper panel shows the pion mass dependence of 

and C14 in the isovector (upper panel) and isosinglet (lower the lattice vector meson (p) mass. The lower panel shows the 

panel) channel, cf. Eq. (17 1, to the lattice vector meson (p) dipole mass m D , cf. Eq. |i"6| , over the lattice vector meson 

mass as, a function of ra^. (p) mass as a function of rn w , in the isovector channel. 



data. Clearly, such a detailed analysis of the Q 2 -dependence of F\ and F 2 is not possible on the basis of the currently 



available lattice data. Still, to get some first insight into the physics behind our preferred parametrization in Eq. (17 1, 
we now explore a matching to a simplistic vector meson exchange ansatz of the generic form 



F{Q 2 



a 1 



a 2 



Ml 



Mf 



Q 2 



(20) 



In the isosinglet channel, one might expect that the lower of the two masses, say Mi, corresponds to the w(782). 
In contrast, in the isovector channel the two-pion continuum contribution plays a leading role, which also generates 
a p(770)-meson exchange contribution. A comparison with the simple ansatz Eq. (20 1 might therefore show that 
M?- d ~ TOp. 



To facilitate the matching of the two-parameter ansatz in Eq. (17) with Eq. (20 1, we implement, in addition to 



a u ~' 
3 a j 



charge conservation, i.e. ^ 
obtained from perturbative QCD 



/(Mf- d ) 2 = 1 for F[ l - d and similar for the isosinglet case, also the large-Q 2 -behavior 



i.e. 



F t {Q 2 ) 



1 / (Q 2 ) l+1 , by setting J^. aj = for the Dirac form factor. 



In the case of F 2 , we would have the additional condition J2j a jMj = 0. We then compute the lowest real solution 
for M\ from the parameters C12 and C14, which were obtai ned from the fits discussed above. The numerical values for 
the extracted masses are provided in Appendix [a] In Fig. 14 we display them in the form of ratios M^~ d /my and 



M 1 l+a /my with m v -- 



(m w ) as functions of the pion mass. As we do not have results available for m£ , we use 

also holds at larger pion masses. Remarkably, 



J at 



,lat 



mjf* instead also in the isosinglet channel, expecting that 
we find that the ratios are very close to, and in most cases within errors fully compatible with, unity over the full 
range of pion masses from m n ~ 1.5 GeV down to m n ~ 0.25 GeV. That this is a non-trivial observation is supported 
by the fact that the lattice vector meson (p) mass, which has been obtained independently, shows a strong pion mass 
dependence, as illustrated in the upper panel of Fig. 15 This pion mass dependence is clearly compensated to a good 



approximation in the ratios in Fig. 14 At the same time, we see from the lower panel in Fig. 15 that the pion mass 
dependences do not cancel out in the ratio mu/mv of the dipole mass, obtained from dipole fits to F\, Eq. (16), to 
the lattice vector meson mass. 

We interpret these results as providing strong evidence for the assumption that the Q 2 -dependence of Fi (within 
the accessible ranges) is to a significant extent driven by vector meson exchange contributions, in particular from 
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FIG. 16. Isovector anomalous magnetic mom ent versus m 



as obtained from fits to F% using Eq. ( 19 1 



FIG. 17. Isosinglet anomalous magnetic moment versus m 7 
as obtained from separate fits to F% and F$ using Eq. ( 19 I 



the u> and p mesons. Furthermore, these findings provide additional support in favor of our preferred 2-parameter 
parametrization in Eq. (17). 



E. Anomalous magnetic moment 



We now turn to a discussion of the anomalous magnetic moment, k — ^(0). As it cannot be extracted directly 
at Q 2 — from a calculation of the Pauli form factor on the lattice with our methods, we have to rely on the Q 2 - 



parametrizations discussed in section IIIB For the reasons given above, we will focus here on the results from the 
more flexible 3-parameter parametrization in Eq. (19). 



Our results for K u -d as a function of the pion mass are displayed in Fig. |16[ While the data points are systematically 
rising as we approach lower pion masses, they are still about 25% below the precisely known experimental value of 
K„,_d = 3.7058893 at the lowest accessible pion masses of ~ 200 — 300 MeV. 

Since the magnitude of is much 



Corresponding results for the isosinglet (u + d) channel are shown in Fig. 



smaller than that of F£ d , the respective lattice data points are in many cases very close to or even compatible with 
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zero, cf. Fig. |4j making a reliable extrapolation in Q 2 very difficult. We therefore have fitted the contributions from 
up and down quarks separately employing the polynomial ansatz in Eq. (19), and subsequently computed n u +d (as 

are again compatible 



17 



well as (f /2 ) from the individual parts. While most of the resulting data points in Fig. 
with zero within uncertainties, we still can observe a systematic trend towards negative values at lower pion masses. 
In the region m n < 500 MeV, we even see an overlap with the experimental value within uncertainties. 

We will take a closer look at the pion mass dependence of k below in section |IVB| As before, we do not find any 
indications for statistically significant systematic discretization or finite volume effects for this observable, as will be 
discussed in section IV CI 



Pauli radius 



The Pauli radius, (r 2 )2, is given by the slope of ^(Q 2 ) at zero momentum transfer. Since the lowest values of Q 2 
for which we can access are in the range of Q 2 ~ 0.15, . . . ,0.5 GeV 2 (for standard periodic boundary conditions 
in spatial directions, and depending on the lattice parameters), the computation of the slope heavily relies on the 
employed parametrization of the Q 2 -dependence. The results for our preferred polynomial, Eq. (19), ansatz are 
displayed in Fig. [18] for the isovector case. Overall, we find that the central values for the polynomial parametrization 
are higher than for the tripole ansatz. At the lowest accessible pion masses, the results from the polynomial ansatz 
in Fig. 18 are just about 20% below the phenomenological value. In contrast, one finds that the corresponding 
lattice data points from the tripole parametrization are about 40 — 50% below the phenomenological number. Not 
surprisingly, the uncertainties from the more flexible 3-parameter fits are significantly larger, and potentially more 
realistic, than for the 2-parameter tripole ansatz. As in section IIIB above, we prefer also in this case the more general 
polynomial ansatz over the standard dipole or tripole parametrizations. With respect to -F^Q 2 ), however, a more 
conclusive assessment probably has to be based on lattice results obtained in larger volumes or employing (partially) 
twisted boundary conditions in order to get access to lower and more densely spaced values of Q 2 . 
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With respect to the isosinglet channel, we first note that (r 2 ) 2 can be written as (r 2 ) 2 — — 6/9 2 /k, where p 2 = 
dF2(Q 2 ) I dQ 2 \Q2 =Q is the slope of the Pauli form factor. Since K u +d turns out to be small and mostly compatible with 
zero within errors over a wide range of pion masses, cf. Fig. [TTJ we will avoid the resulting substantial uncertainties 
in (?" 2 ) 2 +d by considering instead the slope alone, — 6p 2 = (k X {r 2 ) 2 ) u+d . 

It is also interesting to note that for the isosinglet Pauli radius, or more precisely the slope (k x (r 2 ) 2 )" +d , one 
finds a rather widespread range of values from experiment and phenomenology: The super-convergence approach of 
Ref. [8] gives (with a Dirac charge radius of (r 2 ) p E ~ 0.84 fm 2 that is close to the recent measurement by Pohl et 
al. [TS]) (k x (r 2 ) 2 )" +d = 0.04 ± 0.12 fm 2 , where we have obtained the uncertainty from a standard (uncorrelated) 
error propagation. From the same publication [5] the "Recent determinations" from Table I give (k x (r 2 )2) u+d = 
-0.28 ± 0.52 fm 2 , for (r 2 ) p E ~ 0.88 fm 2 that is closer to the PDG value [H]. An even larger negative value can be 
obtained from the parametrization of Ref. [32], (« x {r 2 )2) u+d = —0.66 ± 0.29 fm 2 (taking into account the error 
correlation matrix). In Fig. 19 we show our lattice results for (k x (r 2 ) 2 ) n+d together with the estimated range 
of phenomenological values. While the lattice data points at large pion masses are mostly close to, and within 
uncertainties compatible with, zero, we observe a trend towards non-zero, negative values below m 700 MeV. 
Accordingly, our results at the lowest pion masses, with (k x (r 2 ) 2 )" +d ~ 0.0, . . . , —0.7, are fully compatible with 
the wide range of values from experiment and phenomenology, indicated by the shaded band. We expect the use 
of (partially) twisted boundary conditions to be of great help in order to pin down the parametrization of the (in- 
dependence of F 2 +d at small Q 2 , which should lead in turn to significantly more precise values for (k x (r 2 ) 2 ) u+rf . 



IV. CHIRAL EXTRAPOLATIONS 



In the literature, baryon ChPT has often been applied to lattice results with the goal to extrapolate them from pion 
masses ^ 300 MeV downwards to the chiral limit and to obtain in this way "a priori" or "a posteriori" predictions at 
the physical point. In this process, less known or previously unknown low energy constants (LECs) are treated as free 
fit parameters and are thereby determined from the available lattice data points and their respective uncertainties. 
It is important to keep in mind, however, that in order to provide more than a mere parametrization of the m^- 
dependence of the data, such an approach has to rely on the assumptions that (i) the systematic uncertainties of the 
lattice calculation are reasonably well under control, and (ii) that the particular ChPT-formula to the given order is 
applicable at larger pion masses in the first place. While we are able to study the systematic uncertainties to some 
extent directly on the basis of our extensive sets of lattice ensembles and results, cf. section|Vj the latter assumption is 
gcncrically hard to justify, in particular since for most nucleon observables results from (H)BChPT are only available 
at the 1-loop level, and higher order corrections are difficult to quantify. 

Therefore, in this work, we follow a somewhat different path: 

• Well-known constants like qa of f„ will be fixed as usual to either their physical or their chiral limit values, see 
Table [XXJ. The resulting uncertainties due to a variation of the constants between these values will be studied for 
selected observables. 
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• Low energy constants whose values are at least approximately known will be varied in reasonably wide ranges 
to assess the related uncertainties. 

• Central physical quantities of interest in this work, in particular the chiral limit value of the anomalous magnetic 
moment, k , as well as regularization scale dependent counter term parameters, will be treated as free fit 
parameters. They will be determined by fits preferably only to the experimental and phcnomenological values 
of the observable under consideration at the physical point. Only if this turns out to be insufficient, we will 
include our lattice data points for pion masses below 260 MeV in the fit. 

In essence, we attempt an "upwards extrapolation" from m p h ys towards lattice data points at larger pion masses. In 
combination with our study of potential systematic uncertainties, this approach provides an opportunity to assess 
the applicability of the different available ChPT-schemes in the range in between the physical pion mass and typical 
lowest lattice pion masses of m w ~ 200, . . . , 400 MeV. 

A. Dirac radius 

In [251 [33] , the ttItt and Q 2 -dependence of the nucleon vector form factors was studied in the small scale expansion 
(SSE), a heavy baryon scheme with explicit A-degrees of freedom, to C(e 3 )- The resulting pion mass dependence of 
(r 2 )" -d is given by 

2\«-d,sSE J (- - i /- ~ ■> - /rrw\ > \ c 2 [ — fm 



^ /"x = - — — 1 + 7g A + (l0g A + 2 In ^ + — f- 26 + 30 In 




A 



\J Sm 2 — m 2 



(21) 



which depends on four LECs, the pion decay constant /„■, the isovector axial vector coupling constant g&, the axial 
vector pion- nucleon- A coupling constant ca = g^NA, and the A-nucleon mass difference dm = -m^, as well as a 
counter term b[ t J (A) that removes the regularization scale dependenc^J Generically, the LECs are taken in the chiral 
l imi t, i.e. /„• = /° etc., however to the order considered, they can as well be taken at the physical point. Equation 

l— d 



(21 1 shows explicitly the well-known logarithmic lnm^ divergence that is expected in the chiral limit of 



From Eq. (21), the leading 1-loop HBChPT result (see, e.g., [IHIEQ]) can be easily recovered by setting ca — 0, 
giving 

<'- 3 >r iHB ° hPT - - {1 + TA + (lOd + 2) 1» (^) } + - (22) 



We have employed both the SSE result in Eq. (21 1 as well as the HBChPT expression in Eq. (22 1 to extrapolate from 



(r) 

the physical pion mass upwards in towards the lattice data points. The counter term B w was in both cases fitted 
to the average phenomcnological value at the physical point. This was done for a range of values of the low energy 
constants f v , gA, tun and Sm, which have been varied in between their physical and chiral limit values, cf. Table |TT] 
The coupling ca in the SSE approach has been varied at the same time in the range of ca = 1, . . . , 1.5. The outcome 
of this procedure is shown in Fig. |20[ where the dashed lines outline the uncertainty band from the heavy baryon 
fits, and the shaded band represents the SSE approach. Here and below, a lighter shading is used for the ChPT- 
extrapolation band for pion masses larger than those included in the fit. It is interesting to see that both approaches 
show a rapidly decreasing isovector Dirac radius as the pion mass increases, even leading to an overlap with the lattice 

(r) 

data points at « 250, . . . , 300 MeV. We find that the adjusted counter term parameter B{ ' varies significantly 

(r) 

for the different combinations of parameters and ChPT-approaches: In the SSE-approach, B[ Q sa —1.28, . . . , —0.37, 
while for the HBChPT case, B[ r ] « -0.16, ...,0.0, for a regularization scale of A = 0.89 GeV. As has already been 
noted in [53], it seems doubtful that these two ChPT approaches to the given orders are quantitatively applicable at 
or above the physical pion mass. The overlap with the lattice data points should therefore be interpreted with some 
care, as it might be accidental and not the result of a physically meaningful chiral extrapolation. 



2 Here and below, an analytic continuation of the form (r 2 — l)^ 1 / 2 ^ (j 1 + \/r 2 — 1^ — > =p(l — r 2 )^ 1 / 2 arccos(r) with r = Sm/m^ 
regarded as implicit for > 8m, i.e. r < 1. 
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[ MeV] 


U [MeV] 


9A 


m N [ MeV] 


5m [ MeV] 





86 


1.2 


890 


330 


139 


92 


1.269 


938 


271 



TABLE II. Standard low energy constants at the physical point and in the chiral limit (estimated). We denote the A-nucleon 
mass difference by 5m = mA — mjv . 
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FIG. 20. Pion mass dependence of the isovec tor Dirac 
radius, as obtained from fits to F™~' 
band outlined by the dashed curves, 
represent heavy baryon and SSE 
respectively. For the details, see Eqs. 
surrounding text. 
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chiral extrapolations, 
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21. Pion mass dependence of the isovector Di rac ra- 



as obtained from fits to F 1 using Eq. (17 1. The 



FIG. 

dius, 

shaded band represents a BChPT-fit to an average of the phe- 
nomenological and experimental values at the physical point 
with the counter term as the only fit parameter. The dotted 
line represents the heavy-baryon limit of the central covariant 
fit result. For the details, see Eq. (231 and the surrounding 
text. 



Results for the pion mass dependence of 



A-DOFs) have been obtained in [44]. To <D(p 4 ), it reads 



in the covariant BChPT-scheme of Ref. [35] (without explicit 



2\u-d, BChPT 



Bci + (rt) 



2\u- 



-4(3) + ( r 2)«- 



d,(4) 



(23) 



where the individual higher order contributions (r\) u d <( 3 > 4 ) are given in Appendix |b] in Eqs. (B2| to (B4|. It is 
interesting to note that in contrast to the SSE-expansion to 0(e 3 ), the C(p 4 )-contribution in Eq. (B4| introduces 



a dependence on the coupling cq, which determines the isovector anomalous mag neti c m ome nt in the chiral limit, 
K°' u ~ d =C(i, as we will see below. We also note that the regularization scale in Eqs. (B2) to (B4) has been set equal to 
the nucleon mass in the chiral limit, m° N ~ 0.89 GeV. 

For the covariant BChPT extrapolation, we have varied the LEC cq in a range of 4, . . . , 6 to account for the related 
systematic uncertainties. As before, for a given value of Cq, we have fitted the counter term B cl to the average of 
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the phenomenological values at mP hys , giving a relatively stable B c i(r 
by the error band in Fig. 
Compared to Fig. 
at m n ~ 260 MeV" 



20 



N J 



-1.5 GcV~ . The result is represented 
which in addition includes a variation of the standard LECs as described above, 
the extrapolation band falls off more slowly, and lies about 20% above the lattice data points 
Notably, a heavy baryon expansion of the covariant BChPT result leads to a curve that quickly 
bends upwards above the physical pion mass, as illustrated by the dotted line in Fig. [21] Assuming that the BChPT 
formula to the given order is applicable at the physical pion mass, this would indicate that the corresponding HBChPT 
expansion has a much smaller radius of convergence, and starts to break down already above ~ 100 MeV. In 
comparing HBChPT, SSE and BChPT expansions it is interesting to note that progress with respect to the isovector 
nucleon form factors has been reported very recently in BChPT including the A-resonance in the so-called (5-power 
counting scheme, see [45] and references therein. 

At the one-loop level in HBChPT, with or without explicit A-DOFs, the isosinglet Dirac radius turns out to be 



independent of m v , (r 2 )" 



2\u+d, HBChPT, SSE 



= const. 20, 21 . Furthermore, while the analytical expression of the BChPT 



result of Ref. [U] for (r 2 )" +d shows at first sight a rather non-trivial pion mass dependence, it turns out to be nearly 
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FIG. 22. SSE chiral extrapolation of the isovector anomalous 
magnetic moment. Th e shaded error band represents the 
fit of the SSE Eq. (24 1 with two free parameters (k° and a 



counter-term) to the experimental value and the lattice data 
for m w < 260 MeV. For the details, see Eq. (24) and the 
surrounding text. 
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FIG. 23. BChPT extrapolation of the isovector anomalous 
mag neti c moment. The shaded error band represents fits of 
Eq. (25 1 with two free parameters (c6=K° and a counter- 
term) to the experimental value and the lattice data for < 
260 MeV. The heavy-baryon limit of the central covariant fit 
result is indicated by the dotted line. For the details, see 



Eq. (25 I and the surrounding text. 



flat in practical To the contrary, a rather strong pion mass dependence for this observable is observed on the lattice 
down to ni-ff ~ 230 MeV, cf. Fig. [TT] Lacking any reason to assume that the pion mass dependence suddenly flattens 
off at the physical pion mass, we conclude that the available ChPT results for (r 2 )" +d are most likely not even 
qualitatively applicable at or above mP hys . We therefore refrain from extracting the relevant LECs from fits to the 
phenomenological values at the physical point. It is interesting to note that, although lacking a theoretical foundation, 
a naive linear extrapolation in of the lattice data points below m w = 500 MeV would get reasonably close to the 



experimental and phenomenological values at the physical point in Fig. 11 



B. Anomalous magnetic moment 



The pion mass dependence of K u -d in the small scale expansion to C(e 3 ) can be written as 



„SSE 



i—d 



K u _ d (m 7r ) - 8E[ r \\)m N ml . 



(24) 



where we provide the explicit expression for ~K u -d(rn. n ) in Appendix[Bnn Eq. (Bl ). In addition to the LECs that were 



already discussed above, depends on the isovector anomalous magnetic moment K? u _ d and the isovector nucleon- 

A coupling constant cy = c v in the chiral limit. The counter-term parameter E± (A) removes the regularization scale 
dependence to the given order. Neither cv nor ca are known to great precision. In order to reduce the number of fit 
parameters, we will keep them fixed but perform various fits for cy = —1.5, . . . , —3.5 GeV -1 and ca — 1-0, . . . , 1.5 
in order to assess the related uncertainties. The two remaining unknowns, and e[ (A), are treated as free 

fit parameters and can be obtained from fits to the experimental value at the physical point and the lattice data 
point at ~ 260 MeV. We show the results of the SSE-fits to K u -d in Fig. 22 In the chiral limit, we obtain 



5.3, 



,5.5, which is remarkable 40% above the precisely known value at the physical point. For the counter- 
term parameter, we find values of ~ —2.9, . . . , —5.2 GeV~ 3 for A = 0.89 GeV. Since the extrapolation band 
continues to fall off above m n ~ 300 MeV, it misses the lattice data points to the right, which show little pion mass 
dependence between m T ~ 300 MeV and m w ~ 700 MeV. 

The corresponding expression in the covariant BChPT-approach of Ref. |Hj reads 



.BChPT 



An) 



m 



N 



l6m%mle r we (X) + (K u -d) (3) + (n u -d) 



(4) 



(25) 



Assuming that the anomalous magnetic moment in the chiral limit fulfills > w K u+d ~ —0.36, which is confirmed by the observed 

pion mass dependence of and the corresponding chiral extrapolations, see, e.g., Ref. I26| and our discussion below in section [IV B| 
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FIG. 24. SSE chiral extrapolation of the isosinglet anomalous 
magnetic moment. Th e shaded error band represents the 
fit of the SSE Eq. (26 1 with two free parameters (k° and a 



counter-term) to the experimental value and the lattice data 
for < 260 MeV. 



FIG. 25. BChPT extrapolation of the isosinglet anomalous 
mag neti c moment. The shaded error band represents fits of 
Eq. ( 27 1 with two free parameters (Ks=K^ +d /3 and a counter- 
term) to the experimental value and the lattice data for < 
260 MeV. The dotted line illustrates the heavy-baryon limit 
of the central covariant fit result. 



m 



(n) 



JV 



m 



(») 

JV 
phys 



where m denotes the nucleon mass used in front of F2 in the parametrization of the current in Eq. 



0.938 GeV, and we will explicitly replace by m^ ys in the following. The contributions at 0(p 3 ) 
and 0(p 4 ), i.e. (K u -dp 3 ' 4 ' , are given in Appendix [bJ in Eqs. (B5) and (B6), respectively. 



3b. In our case, 



JV 



These expressions depend on Cq=k° as well as the additional LEC C4, plus a counter-term parameter e\ 06 . Varying 
C4 in the range of 3.2, . . . , 4.0 GeV -1 , we have determined cq and the counter-term from fits to the experimental value 
and the lattice data at it, 



260 MeV. The result is shown in Fig. 23 



in the chiral limit, 



i—d 



Again, we find a rather large value for K u -d 
,5.1, somewhat below the values of the SSE extrapolation. For the counter-term 

3 



22 



for 



1.0 GeV" J for A = 0.89 GeV. The error band is close to the one in Fig. 
' 300 MeV, but then falls off more strongly, already lying a factor of about two 
More interesting is the observation that the non-relativistic limit, ton — > 00, of 



parameter, we obtain e r lm ~ 0.5, . . . 
the SSE case up to pion masses of < 
below the data point at ~ 400 MeV 
the covariant fit, indicated by the dotted line in Fig. 23 drops off even more strongly and starts to deviate from the 
full result already at the physical point. In the case that the covariant approach is at all quantitatively applicable in 
these ranges of the pion mass, this would suggest in turn that the range of applicability of the corresponding heavy 
baryon expansion is much more limited. 

Turning our attention to the isosinglet channel, we first note that in HBChPT at one-loop level, a pion mass 
dependence is only observed in the modified SSE counting scheme of Ref. [22] (denoted by 'scheme C'). It is given by 



K SSE 
K u+d 



u+d - 2AE 2 m N mi , (26) 

where E2 is a counter-term parameter. The result of a fit to the experimental value at the physical point and the 
lattice data point at ~ 260 MeV is illustrated by the shaded error band in Fig. 24 We find a non-zero, negative 
value for the isosinglet anomalous magnetic moment in the chiral limit, /■ 



= —0.40 ± 0.05, while the counter-term 



parameter turns out to be small and compatible with zero within errors, £2 = 
The pion mass dependence in the BChPT calculation of Ref. 44] reads 



-0.08 ±0.12 GeV" 



„BChPT 
l u+d 



phys 



(27) 



where the contributions (K„+d)' 3 ' 4 ' of 0(p z ) and 0(p 4 ) are given in Eqs. (B7) and (B8| in Appendix |b| One finds 
that n u +d depends on the chiral limit value n° u+d = 3k® = k® + , and a counter-term parameter ei 05 (Aj7 As before, 
the regularization scale has been set to A = 0.89 GeV. From a fit to the experimental value and the lattice data point 



at to 
of e 



260 MeV, we obtain a negative 



"u+d 



10.1 



-0.6, . . . , —0.49, with a small value for the counter-term parameter 

With just two 
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0.58, . . . ,0.79 GeV a . The result of the fit is shown by the shaded error band in Fig. 
data points constraining the fit, the band quickly broadens at larger pion masses, thereby prohibiting a quantitative 
assessment. We note, however, that the center of the band provides a good description of the lattice data points up to 
~ 500 MeV. In the heavy-baryon limit, we find a strongly downwards bending curve directly above the physical 
pion mass, illustrated by the dotted line. This indicates once more that the radius of convergence of the heavy-baryon 
approach at one-loop level is limited to the region below m,P hys (assuming that the BChPT result is applicable up to 
TOtt ~ 260 MeV in our fit in the first place). 
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C. Pauli radius 

To separate the pion mass dependence of the slope of ^(Q 2 ) from that of F^O) = k, we focus here on the product 
k x (r 2 )2 instead of (r 2 )2 (see also [26])- This avoids in particular a potential issue related to the expression for «(m ff ) 
that is used in the denominator of the chiral expansion of (r 2 ) 2 = — 6p2(w 7r )/K(m 7r ), where pi = dF2{Q 2 ) / dQ 2 \q2 =Q 
is the slope. Depending on the order of the ChPT-calculation, one is in general allowed to employ ChPT-expressions 
of different orders for fc(m OT ) in the denominator without affecting the overall consistency of the chiral expansion of 
(r 2 ) 2 to the given order. Since the pion mass dependence of k can be rather strong (as discussed in the previous 
section, where K„_j(m>J tys ) increases by as much as 40% as — > 0), the ambiguity in the choice of K^m^) could 
have a significant impact on the uncertainty of the chiral extrapolation of (r 2 )2- 

We therefore consider in the following the SSE and BChPT expansions of K x (r 2 ) 2 . The result in the SSE [22|, 123] 
to C(e 3 ) is 

U x (, 2 ) 2 r d ' SSE = + c lmN ln( 6rn + J** - l) + 24m N B c2 , (28) 

V ) 8nf%m« 9TT 2 f 2 ^/Sm 2 -ml \m n y m 2 J 

showing explicitly the well-known linear divergence in m w expected in the chiral limit. It depends again on the 
coupling ca, and a counter-term parameter B C 2- As in the case of K u -d, we have varied ca — 1.0, . . . , 1.5, and 
determined B C 2 as the only free parameter from a fit to the phenomenological value at the physical point. We find 



B c2 — 0.32, . . . , 0.82 GeV 3 . The result is illustrated in Fig. 126] by the shaded bands, in comparison to the lattice 
data points obtained from the polynomial ansatz Eq. (19) for the Q 2 -dependence of F 2 . While the extrapolation 



band quickly decreases above TOP hys , it overshoots the lattice results by about 20 — 40% in the region of m„ 
260, . . . , 500 MeV. 

A significantly more involved expression for the m^-dependence has been obtained in the BChPT-scheme of Rcf . 
to 0(p 4 ), which can be written as [H] 



(« x (r 2 ) 2 ) = (24m° N e r 74 (X) + {nr 2 )^) + (nr 2 )-^) , (29) 



where the individual terms, {nr 2 ) u and (/trf)" "U 4 ), are provided in Appendix b] in Eqs. (B9) and (B10 1 



This result depends, apart from the counter-term parameter 674(A), also on the couplings c 4 and C6=K^_ d . The 
regularization scale has been fixed to A = 0.89 GeV. To study the predicted m^-dependence, we have varied, as 
before, c 4 = 3.2, ... ,4.0 GeV -1 , and used the corresponding values obtained for c 6 = from the BChPT analysis 
of above. The unknown parameter e^ 4 has been fitted to the phenomenological value at the physical point, 

giving ey 4 ~ 1.6,..., 2.0 GeV -3 . We com par e this approach with the m^-dependence of our results obtained from the 
polynomial ansatz for F2~ d {Q 2 ) in Fig. I27J Similarly to the SSE-extrapolation discussed before, the extrapolation 
curves lie somewhat above the lattice data points for ~ 260 MeV to 500 MeV. Taking the heavy-baryon limit of 
the central band, we obtain the dotted curve in Fig. [27] Under the assumption that the BChPT result is applicable 
at the physical pion mass, we find that the contributions of C(l/(mjv)™), which are included in the covariant BChPT 
approach, start to play an important role already before mP hys is reached. Hence also for this observable, our results 
indicate that the range of applicability of the leading order heavy-baryon ChPT result is restricted to pion masses 
< 100 MeV. 

With respect to the isosinglet channel, we note that at the one- loop level in HBChPT, (r 2 )'^ +d is predicted to 
vanish, i.e. the form factor is independent of Q 2 , F 2 li+d (Q 2 ) = const. [5T]. Similar to the case of the isosinglet Dirac 
radius, the BChPT calculation gives at first sight a rather non-trivial "m-Tr-dependence [44] . but in practice it turns out 
to be flat over the full range of relevant pion masses, (n x ( r 2 } 2 ) tl + ,:i . BChPT ~ const. Given the poorly determined value 
at the physical point, we conclude that it is currently difficult to provide even a semi-quantitative chiral extrapolation 
of (k x (r 2 )2) u+d from mP hys to the chiral limit and to larger pion masses. 

V. SYSTEMATIC UNCERTAINTIES 

A. Contaminations from excited states 

A potentially important source of systematic uncertainties is given by contributions from excited states in the 
nucleon correlation functions. For too small distances between the operator insertion time r and the source and sink 
times they could adversely affect the ratios of three- to two-point functions in Eq. Q , and thereby the plateau values 
from which we extract the form factors. In turn, if the sink time i sn k and accordingly r are chosen too large, the 
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FIG. 26. SSE chiral extrapolation of k x (r 2 }2 in the isovector 
channel. The shaded error bands represent fits of the SSE 
results in Eq. ( |28| l with one free parameter (a counter-term) 
to the experimental value. 



FIG. 27. BChPT extrapolation ofcx (r 2 )2 in the isovector 
channel. The shaded error bands represent fits of the BChPT 
results in Eq. ( |29[ ) with one free parameter (a counter-term) 
to the experimental value. 



signal-to-noise ratio begins to deteriorate, and the data points start to fluctuate more strongljj^] This is mostly an 
issue at larger hadron momenta required for the analysis of the form factors at Q 2 > 0. Hence we have to seek a 
compromise between potential contaminations from excited states on the one hand, and noisy /fluctuating correlation 
functions and plateaus on the other. 

In this work, we have chosen primarily a fixed distance between source and sink of about 0.95 fm. Instead of studying 
the excited states contributions directly by performing, e.g., multi-exponential fits (which are notoriously unstable) 
of the correlations functions, we have analyzed the form factors for a range of different sink times i S nk = 11, . . . , 19, 
for a single ensemble with j3 — 5.29, n — 0.13590. The dependence of, e.g., the data points for F\ at fixed Q 2 on i sn k 
is then a direct indicator for the possible influence of excited state contributions on our results. The results of this 



study are displayed in Fig. 28 showing F" _d as a function of t sn k for the four values Q 2 ~ 0.49, 0.94, 1.36, 1.75 GeV 2 . 
We note that the broader band (corresponding to our primary choice i sn k = 13) is compatible with all data points 
at sink times up to and including i sn k = 16 within errors. While the central values decrease on average by a small 
amount as i sn k increases from 13 to 16, no clear systematic trend can be established when the uncertainties are taken 
into account. At large i sn k > 17, we find that the data points start to fluctuate more strongly as the momentum 
transfer increases. This indicates that the plateaus indeed become unstable due to deteriorating signal-to-noise ratios 
of the correlation functions in the ratio at large times and momenta. At larger Q 2 > 1.8 GeV 2 (not shown), we even 
find that the extracted values for F\ quickly approach zero as t sn k — > 19. 

In summary, for the given ensemble, we cannot identify a systematic dependence of our results for F\ on the 
sink time within errors, excluding large i sn k values where strong fluctuations and low signal-to-noise ratios make a 
quantitative analysis impossible. This indicates that the uncertainty due to excited state contaminations is not larger 
than the statistical errors in our study. 

B. Discretization effects and pion mass dependence at fixed Q 2 

Studies of discretization effects and the pion mass dependence of the lattice results are usually directly performed 
for the fundamental observables of interest. In our case, these are the radii of the Dirac and Pauli form factor, as 
well as the anomalous magnetic moments. However, due to the discrete values of the momentum transfer that can 
be accessed in a finite volume, in particular the still rather large, lowest non-zero Q 2 of « 0.2 GeV 2 in our case, the 
extraction of these observables requires non-trivial inter- and extrapolations of the form factor data in Q 2 . To avoid 
an intermixture of the primary lattice artifacts with uncertainties due to the required Q 2 -parametrization, we now 
attempt to investigate the a- and m^-dependences of the form factor data directly for fixed values of Q 2 > 0. In this 
regard, one has to keep in mind that changes in the lattice volume, spacing, and the quark mass (and thereby the 
nucleon mass), lead in general to different sets of values of Q 2 for the different ensembles. To study the lattice spacing 
dependence, we have therefore scanned our data sets for narrow ranges in Q 2 and (with maximum relative widths 



4 Also the two-point functions taken at t sn k in the ratio might fluctuate around zero within errors, thereby leading to an unreliable result 
for the plateau. 
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FIG. 28. Dependence of F™~ at fixed values of Q 2 on the sink time i S nk in the three point function, for /3 = 5.29, k = 0.13590. 
Results for our choice t sn k = 13 are indicated by the filled (red) points and corresponding bands. The thin gray shaded band 
represents the parametrization of Ref. |H] of the experimental data. 
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FIG. 29. Lattice spacing dependence of the isovector Dirac form factor for fixed narrow ranges in Q 2 and m w . The gray shaded 
bands represent the parametrization of Ref. [32] of the experimental data. 



of 8%) for which data points for three or more couplings j3 are available. The residual Q 2 - and m^-dependences within 
these narrow windows were taken into account by interpolations and subsequent relative shifts of the data points to 
the central values of Q 2 and in the respective ranges. As a test, we have monitored the relative shifts and found 
that their absolute values are about the same size as the statistical errors of the shifted data points, 
independence are displayed in Figs. 



The results for the 
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fluctuations of the central values as functions of a arc visible, t 



and 



30 



for F%- 



and ir- 



respectively. While some 



icy do not seem to follow a systematic pattern. Overall, 



the data points are compatible with a constant behavior within statistical errors. In combination, the uncertainties 
and the fluctuations of the data for the given Q 2 -values are however too large to allow for a consistent, quantitative 
continuum extrapolation. Still, although we cannot exclude the presence of some discretization effects, we do not see 
any evidence that they could significantly reduce the large gap between the lattice data points and the experimental 
result illustrated by the shaded bands in Figs. [29] and [30] 

A similar approach to the pion mass dependence at fixed Q 2 leads to the results displayed in Figs. 31 and 32 
Following the above findings on the a 2 -dcpcndcnce, we have in this case included all /3 on an equal footing. For 
F?~ d in Fig. ^ 
up to ~ 1.8 
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we observe an approximately linear dependence on m w over a wide range of Q 2 from ~ 0.5 GeV z 
GeV 2 . While the data points do show a slight downwards trend in the right direction, simple linear 
extrapolations would clearly miss the experimental values by about 20% to 40% at the physical pion mass. Keeping 
in mind that chiral perturbation theory predicts a logarithmically diverging slope of F^~ d at Q 2 — as —> 0, it 
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FIG. 30. Lattice spacing dependence of the isovector Pauli form factor for fixed narrow ranges in Q 2 and m^. The gray shaded 
bands represent the parametrization of Ref. 32 of the experimental data. 
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FIG. 31. Pion mass dependence of the isovector Dirac form factor for selected ranges in Q 2 and m v . The gray bars represent 
the parametrization of Ref. [32] of the experimental data at the physical point. 



is not surprising that also for d (Q 2 ^ 0), a non-linear m^-dependence has to set in at low pion masses. 



In the case of F% d in Fig. 32 we find again that the lattice results are, to a good approximation, linear in m ff . 



While the data points at lower Q 2 are constant within the uncertainties, a slight downwards slope seems to develop 
as we approach larger momentum transfers. In contrast to F^~ d , a naive linear extrapolation in would even lead 
to an overlap with the experimental values at the physical pion mass, at least for the lowest value of Q 2 ~ 0.50 GeV 2 
in Fig. 32 This does not imply, however, that the lattice results for m T > mP hys provide a good description of the 



experimental data over a wider range of Q 2 . At large Q 2 , linear extrapolations in m n would lead to values for 
that are systematically larger than in experiment. We therefore find again that the Q 2 -slope of the lattice data (even 
when naively extrapolated to m|J hyB at fixed Q 2 ) is too small, and a typical dipole or tripole extrapolation (see section 
III B I to Q 2 = would then lead to a F 2 " _<i:lat ((5 2 = 0) = K^_ d < . The apparently good agreement of the 



lattice data points with the experimental values in the top row of Fig. [32] has to be interpreted as the result of a too 
small slope and, at the same time, a too low normalization (at Q 2 — 0) at unphysically large lattice pion masses. 
This is studied explicitly in sections |IVB| and |IV C| above, where we discuss the pion mass dependence and chiral 
extrapolations of K u -d and the slope of F^ _d (Q^), respectively. 
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FIG. 32. Pion mass dependence of the isovector Pauli form factor for selected ranges in Q and m n . The gray bars represent 
the parametrization of Ref. [35] of the experimental data at the physical point. 
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FIG. 33. Lattice spacing dependence of the isovector Dirac radius for fixed ranges of m n , as obtained from fits to using 
Eq. ([IT]) (cf. Fig. |D . The residual pion mass dependence in the given m^-windows has been taken into account by restricted 
linear fits, and the data points have been shifted accordingly to the central m^-values. The gray bands represent the range of 
values obtained from experiment and phenomenology at the physical point. 



C. Discretization and finite volume effects in the Dirac radius and the anomalous magnetic moment 



Our study of the a-dependence of F\ at fixed m n and Q 2 above in section VB already indicated that the dis- 
cretization errors are small. Here, we perform a si milar analysis for the isovector Dirac radius obtained from the 



we show our results for 



for narrow ranges of as 



Q 2 -parametrization based on Eq. (|17|). In Fig. 

function of a 2 . For direct comparison, we also show an average of the experimental results at the physical point as a 
gray error band. The residual pion mass dependence of the lattice data within the m^-ranges has been accounted for 
by linear fits to the pion mass dependence and subsequent relative shifts of the data points to the central m^-values. 
As expected, we find that the results are compatible within statistical uncertainties for the three or four different 
available values of a 2 . Apart from some small fluctuations, which, however, do not show a systematic trend, we find 
that even the central values of the data points are in good agreement. Overall, in the accessible parameter ranges, 
and for a 2 - 0.0035, . . . , 0.007 fm 2 , we therefore do not observe any significant, systematic lattice spacing dependence 
of our results. Most importantly, a naive extrapolation of our data in a 2 to the continuum limit would not bring us 
any closer to the experimental value indicated by the gray band. 
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FIG. 34. Volume dependence of the isovector Dirac radius for fixed values/ranges of m^, as obtained from fits to F'^~ d using 
Eq. (17 1 (cf. Fig. [8j. The residual pion mass dependence in the given m„- windows has been taken into account by restricted 
linear fits, and the data points have been shifted accordingly to the central m^-values. Fits to the //-dependence of the data 
points are indicated by the dashed lines, and the light shaded bands represent the resulting value in the infinite volume limit. 
The upper gray bands indicate the range of values obtained from experiment and phenomenology at the physical point. 



With respect to finite volume effects, we display in Fig. 34 our results for (r 2 )i as a function of the box length L 



As before, the lattice data points were shifted to the central values of the indicated narrow ranges in m v employing 
linear interpolations in order to account for the residual pion mass dependence. In contrast to the absence of any 
a 2 -dependence discussed before, we observe a slight, systematic upwards trend of the data points as L increases. In 
a first attempt to quantify this observation, we have fitted the L-dependence of the data points in each pion mass 
range with a simple exponential ansatz inspired by predictions from chiral perturbation theory: a + 6exp(— m^L). 
The results of the fits are indicated by the dashed lines, and the corresponding estimated values in the infinite volume 
limit are shown as light shaded bands. We note that all data points, apart from the ones with x L < 3.4 (filled 
diamonds), show at least a small overlap with the infinite volume band within uncertainties. The rightmost points 
at larger volumes are in all cases fully compatible with the estimated results at L = oo and hence can be regarded as 
corresponding to the infinite volume limit. 

We now turn to systematic uncertainties in the anomalous magnetic moment, following the sam e st rategy as outlined 



above for the case of the isovector Dirac radius. The a 2 -dependence of K u -d is shown in Fig. 35 for two ranges of 
m % . Within the uncertainties, the data points do not show any systematic trend as the lattice spacing decreases 
and are fully compatible with constants in a 2 . As before, a linear extrapolation would not bring us any closer to the 
experimental value indicated by the thin gray band. 

Figure 36 displays the dependence of n u -d on the box length L for selected ranges of m^. In contrast to " "' 



in Fig. |34| the data points do not show any clear upward or downward trend as L — >■ oo. Since the uncertainties and 

fluctuations are somewhat larger, we will have to leave a more quantitative estimate of finite volume effects in k u ^ 

for future works. 



VI. SUMMARY AND OUTLOOK 



On the basis of an extensive set of ensembles for two flavors of C(o)-improved Wilson fermions and Wilson gluons, 
we have computed and studied the Dirac, Fi(Q 2 ), and Pauli, ^(Q 2 ), form factor of the nucleon. Four different 
lattice spacings from a = 0.083 fm to a = 0.060 fm, spatial volumes of V ~ (1.0, . . . ,3.0 fm) 3 , and a wide range 
of pion masses extending down to ~ 180 MeV, allowed us to investigate in some detail the continuum, infinite 
volume, and chiral limit. As in previous studies, we do not yet see an overlap or agreement of the lattice data with 
results from experiment and phenomenology for the slopes of the isovector Dirac and Pauli form factors, nor for 
the normalization of the latter, i.e. the isovector anomalous magnetic moment K u -d = F% (Q 2 — 0), even at the 
lowest accessible pion masses of 180 MeV < < 300 MeV. Our results indicate that these discrepancies cannot be 
explained by discretization or finite volume effects. Contributions from excited states, another source of systematic 
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linear fits, and the data points have been shifted accordingly to the central m 
from experiment at the physical point. 



35. Lattice spacing dependence of the isovector anomalous magnetic moment for fixed ranges of m^, as obtained from fits 
~ d using Eq. ( 171. The residual pion mass dependence in the given m^-windows has been taken into account by restricted 

values. 



The gray bands represent the value 





FIG. 36. Volum e de pendence of the isovector anomalous magnetic moment for fixed ranges of m w , as obtained from fits to 
F2~ d using Eq. |l7| . The residual pion mass dependence in the given m^-windows has been taken into account by restricted 
linear fits, and the data points have been shifted accordingly to the central m^-values. The gray bands represent the value 
from experiment at the physical point. 



uncertainties, were studied for a single ensemble by varying the sink time of the three point function, also giving 
no hint for substantial effects larger than the statistical errors. Concerning the pion mass dependence, our results 
for Fi(Q 2 ) and F 2 (Q 2 ) at fixed values of Q 2 look mostly linear in or m 2 within the accessible ranges. Linear 
extrapolations to the physical pion mass, however, would not lead to an agreement with experiment. 

Our data for the ratio {F? /F?)(Q 2 ), for pion masses below ~ 500 MeV, show a reasonable overlap with the 
phenomenological parametrization over the full range of Q 2 we could access. It is interesting to observe that this ratio 
drops off by about 50% reaching Q 2 ~ 2 GeV 2 , pointing towards a much narrower spatial distribution of up quarks 
in the proton than of down quarks. Furthermore, we find that the Dirac radius of down quarks is systematically 
larger than for up quarks, over the full range of available pion masses. Concerning the ratio of the Pauli to the Dirac 
form factor, (Fa/ ))■"' , our results are in general compatible with a rather flat Q 2 -dependence as observed in 
experiment, although the statistics and the covered Q 2_ran g es are a t this point insufficient to permit a quantitative 
assessment. With respect to results in the isosinglet channel or involving individual quark flavors, one has to keep in 
mind that quark line disconnected contributions have been neglected. 

For inter- and extrapolations in the momentum transfer, and in particular to extract the mean square radii and 
anomalous magnetic moments, we have performed and studied different parametrizations of the Q 2 -dependence of F\ 
and -FV For the Dirac form factor, we find that a more flexible (polynomial) 2-parameter ansatz provides a numerically 
and physically much more convincing description compared to the commonly employed dipole fits. This is borne out 
by a matching onto a basic vector meson exchange model: In the case of the polynomial parametrization, we find that 
the extracted lowest vector meson masses agree very well with the separately computed lattice vector meson (p) mass, 
over a very wide range of pion masses from 1500 MeV down to 260 MeV. This indicates that the Q 2 -dependence of 
the nucleon form factors on the lattice is to a significant extent governed by the exchange of the lowest vector meson 
resonances, that is the p in the isovector, and the u in the isosinglet channel. 

With respect to chiral extrapolations using chiral perturbation theory, we followed a somewhat different path than in 
the past. Instead of attempting extrapolations of the lattice data down to the physical pion mass, we investigated the 
applicability of the different ChPT-schemes by including the known results from experiment and phenomenology at 
the physical point, and, only if necessary, lattice data for pion masses below 260 MeV, in the chiral fits. Not precisely 
known low energy constants were varied over sufficiently wide ranges to assess the related uncertainties. While our 
data points show for the first time the onset of a non-analytical chiral behavior at the lowest pion masses, it still turns 
out to be difficult to achieve a consistent quantitative understanding of the m^-dependences using the different heavy 
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baryon, small scale expansion (explicitly including the A resonance), and covariant BChPT approaches. This is in 
particular the case for the isovector Dirac radius, where traditional HBChPT predicts a too steep, and a covariant 
BChPT approach a too flat slope as the pion mass increases above mjj llys . Also the extrapolations of the isovector 
anomalous magnetic moment and slope of the Pauli form factor are still challenging, as they clearly under- or overshoot 
the lattice data points for pion masses of Si 300 MeV. In the chiral limit, we obtain a rather large ~ 4.8, . . . , 5.5, 
in agreement with previous lattice studies. For the anomalous magnetic moment in the isosinglet channel, we obtain 
a reasonable description of the lattice data up to pion masses of ~ 500 MeV, within rather broad ChPT extrapolation 
bands. In this case, we find a clearly negative value of K^ +d ~ —0.6, . . . , —0.35 at m n — 0. Apart from K u +d, we did 
not attempt any chiral extrapolations in the isosinglet channel, since the available 1-loop results from ChPT in this 
sector are clearly not applicable even at the physical pion mass. In all considered cases, we find that the heavy-baryon 
limits of the covariant BChPT extrapolations appear to break down at or even below mP hys . This casts strong doubts 
on the applicability of the leading 1-loop heavy-baryon approaches in this region of the pion mass. 

Recalling the importance of vector meson exchange contributions for the Q 2 -dependence of the form factors that 
we observed before, it could be interesting to include such contributions explicitly in the ChPT-description not only 
of the Q 2 -, but also the m^-dependence of these observables. Form factor calculations including explicit p, uj and <p 
resonances have been performed in covariant BChPT in the so-called EOMS-scheme [4BTI48] . leading to an improved 
description of the Q 2 -dependence of the nucleon form factors at the physical pion mass. It will be interesting to study 
the applicability of such calculations at larger (lattice) pion masses, and eventually to compare with the combined 
TOtt- and Q 2 -dependence of the lattice data. 

Concerning future nucleon form factor studies on the lattice, our current analysis underlines the importance to obtain 
results for pion masses below 200 MeV in sufficiently large volumes of V gl 3.5 fm, which represents a remarkable 
computational challenge. Apart from being crucial for the chiral extrapolation and comparison with experiment, 
such calculations will be indispensable for a quantitative understanding of the volume dependence at our lowest pion 
mass of ~ 180 MeV, where we begin to see finite size effects in, e.g., the data for the Dirac radius. Furthermore, 
our investigation and comparison of different ansatze for the Q 2 ~dependences has shown that precise data points are 
required over a wide range of the momentum transfer in order to limit additional parametrization uncertainties, in 
particular for the extraction of the anomalous magnetic moments and the radii from the slopes at Q 2 = 0. In this 
respect, (partially) twisted boundary conditions for the quark fields in spatial directions have already proven to be 
highly helpful to access very small non-zero values of the momentum transfer in the case of the pion form factor, see, 
e.g., [39H5T] . First studies along these lines for the nucleon form factors are promising [53] and will be continued in 
the near future. Regarding higher Q 2 > 2 GeV 2 (involving larger nucleon momenta), it will be important to carefully 
monitor fluctuations in the correlation functions and potential contaminations from excited state contributions. 



Appendix A: Collection of numerical results 



Table III shows our results for the mean square radii, anomalous magnetic moments, and vect or me son masse s 
M" _d and M'" +d , for all ensembles specified in Table [ij Definitions and details are given in sections 



IIIB 



and 



HID 



Appendix B: Chiral perturbation theory formulae 



Here we provide a collection of (parts of) SSE and covariant BChPT expressions for the mean square radii and 
the anomalous magnetic moments. For the details, we refer to Refs. 22, 2"3"1 133] and the sections IV IV B and IV C 
above. 
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(r 2 >r d [ fm 2 ] 


(r 2 )^ +d [ fm 2 ] 


K u — d 




<r 2 >r d [ fm 2 ] 


(k X (r 2 ) 2 y +d [ fm 2 ] 


M?~ d { GeV] M? +d [ GeV] 


1 


0.103(11) 


0.125(8) 


1.461(148) 


-0.121(155) 


0.145(41) 


-0.147(111) 


1.602(209) 


1.629(133) 


2 


0.186(20) 


0.207(15) 


2.096(909) 




0.320(273) 




1.194(163) 


1.307(102) 


3 


0.194(33) 


0.236(29) 


2.165(2.120) 








1.135(190) 


1.250(102) 


4 


0.100(11) 


0.135(10) 


1.507(273) 


0.406(289) 


0.179(79) 


0.233(262) 


1.542(236) 


1.647(158) 


5 


0.127(14) 


0.144(17) 


1.801(1.562) 








1.485(207) 


0.962(119) 


6 


0.200(6) 


0.267(5) 


2.107(68) 


-0.131(85) 


0.343(23) 


-0.195(90) 


1.129(30) 


1.131(50) 


7 


0.230(10) 


0.322(9) 


2.200(135) 


-0.260(158) 


0.367(45) 


-0.242(179) 


1.020(40) 


0.986(65) 


8 


0.083(6) 


0.111(5) 


1.425(138) 


0.164(119) 


0.158(36) 


0.040(67) 


2.109(205) 


1.885(155) 


9 


0.125(9) 


0.141(9) 


1.676(347) 


-0.498(306) 


0.215(93) 


-0.401(290) 


1.578(174) 


1.428(141) 


10 


0.353(132) 


0.516(173) 










0.789(188) 


0.659(191) 


11 


0.128(15) 


0.170(16) 


1.222(501) 








1.692(193) 


1.252(146) 


12 


0.160(5) 


0.198(3) 


1.924(49) 


0.080(44) 


0.279(17) 


0.028(43) 


1.325(48) 


1.438(40) 


13 






7.885(3.210) 








1.027(223) 


1.017(165) 


14 


0.160(25) 


0.223(15) 










1.257(160) 


1.009(37) 


15 


0.193(4) 


0.259(3) 


2.163(59) 


-0.109(58) 


0.374(21) 


-0.104(63) 


1.146(19) 


1.128(27) 


16 


0.237(10) 


0.326(11) 


2.235(143) 


-0.434(171) 


0.381(50) 


-0.457(209) 


1.010(33) 


0.972(66) 


17 


0.250(21) 


0.445(19) 


2.396(491) 




0.482(148) 




1.064(27) 


1.008(55) 


18 


0.296(20) 


0.417(14) 


2.877(184) 


0.005(273) 


0.602(84) 




0.863(84) 


0.828(86) 


19 


0.319(15) 


0.429(7) 


2.466(80) 


-0.400(91) 


0.474(40) 


-0.413(146) 


0.817(30) 


0.797(24) 


20 


0.330(60) 


0.435(58) 


3.475(2.301) 








1.030(98) 


0.870(136) 


21 


0.112(3) 


0.138(2) 


1.513(40) 


0.039(38) 


0.177(11) 


-0.027(24) 


1.660(99) 


1.713(50) 


22 


0.146(6) 


0.183(5) 


1.792(74) 


-0.074(90) 


0.254(25) 


-0.088(74) 


1.337(58) 


1.364(96) 


23 


0.169(11) 


0.226(9) 


1.789(164) 


-0.182(200) 


0.230(51) 


-0.197(181) 


1.353(137) 


1.397(84) 


24 


0.221(14) 


0.252(10) 


2.031(369) 




0.335(121) 




1.156(124) 


1.105(56) 


25 


0.198(25) 


0.268(16) 


2.098(1.508) 








1.162(179) 


1.009(12) 


26 


0.219(9) 


0.314(7) 


2.248(99) 


-0.278(100) 


0.399(38) 


-0.240(125) 


1.069(42) 


0.977(37) 


27 


0.215(18) 


0.389(24) 


2.505(311) 




0.428(98) 




1.011(42) 


0.982(48) 


28 


0.299(11) 


0.438(8) 


2.608(115) 


-0.269(133) 


0.548(60) 


-0.146(234) 


0.851(24) 


0.766(14) 


3LE III. Results for the mean square radii and anomalous 


magnetic moments, as well as the vector meson masses M" d 



and 

M" +d obtained from a matching to the vector meson exchange ansatz, cf. section HID Entries with very large uncertainties 
have been replaced by ellipses. The ensembles 1,...,28 are sp ecified in Table |I| All results are based on the polynomial 
parametrizations of Fi(Q 2 ) and F2(Q 2 ) discussed in section IIIB 
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In the following expressions, m? N denotes the nucleon mass in the chiral limit, while m^ represents the pion mass 
dependent nucleon mass, miq{m. n ) (44] . 
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